The methods to deal with missing data, implicitly assume a missing data mechanism.

MCAR: the most strict assumption. In practice it is also easiest to deal with MCAR data.

1. Analyze the observed sample only (this will result in unbiased estimates).
2. Use an imputation method to boost the power the the amount of missing data is too large.

MAR: less strict assumption. Most advanced missing data methods assume this mechanism (e.g. multiple imputation, FIML).

• Include variables in study that may explain the missing data, a MAR assumption may become more plausible (as compared to MNAR).
• These auxiliary variables may also help in dealing with the missing data.
• Auxiliary variables: variables related to the probability of missing data or to the variable with missing data.
  • Can be used as predictors in an imputation model or as covariates in a FIML model to improve estimations.

MNAR: least strict assumption.

• MNAR data are also referred to as non-ignorable, because these cannot be ignored without causing bias in results. MNAR data are more challenging to deal with.
• Sensitivity analyses to explore the extend to which the effect is robust against the MNAR data.

• Most important: assumes that the imputed value is an actual observed value.
• Some methods assume MCAR (i.e. regression imputation, random value imputation)

Other examples:

• All missing values are average (i.e. mean imputation)
• People with missing data remain stable (i.e. LOCF)

• Imputing the average value for all missing entries
• Very easy method

Parameter estimates:

• Mean is unbiased when MCAR
• Regression weight/correlation never unbiased
• Standard errors are underestimated

Conditional mean imputation: estimate the imputed value as a predicted value from a regression model: \[Y_{imp} = \beta_0 + \beta_1 * X\]

Parameter estimates:

• Mean is unbiased when MAR
• Regression weight is unbiased when MAR
• Correlation is never unbiased
• Standard errors are underestimated

Ad hoc method for longitudinal data: use the previous observed value to impute the missing values.

• Assumes that people that drop out of the study remain stable.

Parameter estimates:

• Mean is biased
• Regression is biased
• Correlation is biased
• Standard errors are underestimated

The trajectories over time, with LOCF imputations in red.

Below only the cases with missing observations, separated by group.

• Imputing the missing data entries with multiple “plausible” values.
• Takes imputation uncertainty into account
• A method to improve the main analysis results, so having a fully observed data set is not the goal.
  • Unbiased estimates
  • More precision and thus power
• R package: mice

• Imputation phase
• Analysis phase
• Pooling phase

1. Incomplete data
2. Generate multiple copies of the same dataset, but each time differenty imputed values
3. Analyze each imputed dataset
4. Pool results for analyses to final study result

• Imputing the missing data entries with multiple “plausible” values
• Specify the imputation model
  • method to estimate imputed values
  • predictors used in the model
  • can both be specified for each variable with missing data

Per variable with missing data:

• Impute missing value by using the predicted value estimated from observed data.
• Add residual error to simulate sampling error
• Depending on the variable type that needs to be imputed, the imputation method is defined.

Later more on methods

Imputation model = method + predictors

• For each variable, the most relevant predictors can be indicated
• Use variables that are in the analysis model also in the imputation model (compatibility)
• Add auxiliary variables
• Advice to use no more than ~25 predictors (also depends on sample size)

Auxiliary variables: variables related to the probability of missing data or to the variable with missing data.

• One iteration = one cycle through all variables that need to be imputed.
  • The algorithm starts with a random draw from the observed data.
  • Then imputes the incomplete data in a variable-by-variable fashion.
• Typically 5 to 10 iterations before an imputed data set is generated.
• Iterations are repeated to ensure convergence of the algorithm

Iterations for generating the imputations for Solar.R in the first data set:

Iterations for generating the imputations for Solar.R in the second data set:

Iterations for generating the imputations for Solar.R in all five data sets:

• Changes in imputation model may be needed when:
  • Two variables are highly correlated (multi-collinearity)
  • A variable has so many missing values that imputations are unstable

• Theoretically, higher \(m\) is better.
• Usually, use \(m = 5\) for model building and \(m\) ~ % missing data for final analysis.
• In lower \(m\) the between-imputation-variance is enlarged by \(\frac{1}{m}\) before computing the total variance.

• Each imputed data set is analyzed, with the substantive analysis model

• This results in \(m\) sets of results

• Workflow:

  • Imputation
  • Analysis
  • Pooling

> # A tibble: 10 x 6
>    term        estimate std.error statistic   p.value  nobs
>    <chr>          <dbl>     <dbl>     <dbl>     <dbl> <int>
>  1 (Intercept)  24.8       5.71        4.34 0.0000258   153
>  2 Solar.R       0.0996    0.0278      3.58 0.000455    153
>  3 (Intercept)  21.7       5.62        3.85 0.000173    153
>  4 Solar.R       0.105     0.0271      3.88 0.000153    153
>  5 (Intercept)  21.3       5.53        3.86 0.000165    153
>  6 Solar.R       0.103     0.0267      3.87 0.000164    153
>  7 (Intercept)  21.8       5.55        3.93 0.000127    153
>  8 Solar.R       0.109     0.0271      4.02 0.0000930   153
>  9 (Intercept)  20.0       5.55        3.60 0.000433    153
> 10 Solar.R       0.112     0.0271      4.14 0.0000578   153

• Pool the analysis results to obtain final parameter estimates.
• For normally distributed parameters: Rubin’s Rules

• Pooling of point estimates that are normally distributed over the imputed datasets.
• Means, standard deviations, regression estimates, linear predictors, proportions.
• For pooling point estimates, use mean:

\[\hat\theta = \frac{1}{m} \sum^m_{i=1}{\theta_i}\]

• Pooling of variance or standard error around the estimate: combine the within and between imputation variance.

Between variance:

\[\sigma^2_{between} = \frac{\sum^m_{i=1}(\beta_i - \overline\beta)^2}{m-1}\]

Within variance:

\[\sigma^2_{within} = \frac{\sum^m_{i=1}\sigma^2_i}{m}\]

Total variance:

\[\sigma^2_{total} = \sigma^2_{within} + \sigma^2_{between} + \frac{\sigma^2_{between}}{m}\]

• Use a transformation to approximate a normal distribution
• Examples:
  • Correlation: Fisher \(z\)
  • Odds Ratio, Relative Risk, Hazard ratio: log transformation
  • Explained variance: Fisher \(z\) for \(\sqrt{R^2}\)

• riv: Relative increase in variance due to nonresponse
• df: Residual degrees of freedom for hypothesis testing
• lambda: Proportion of total variance due to missingness
• fmi: Fraction of missing information for the parameter estimates

In the imputation phase, imputed values are estimated using an imputation method.

• Regression imputation: norm.predict
• Stochastic regression imputation (regression + residual error): norm.nob
• Stochastic regression with parameter uncertainty: norm
• Predictive mean matching: pmm

• Imputed value: \(Y_{imp} = \hat{\beta}_0 + X_{mis}\hat{\beta}_1 + \epsilon\)

• Parameters \(\hat{\beta}_0\) and \(\hat{\beta}_1\) are estimated from the observed data.

• And updated after the first imputation with the imputed data.

• Imputed value: \(Y_{imp} = \hat{\beta}_0 + X_{mis}\hat{\beta}_1 + \epsilon\)

• Parameters \(\hat{\beta}_0\) and \(\hat{\beta}_1\) are estimated from the observed data.

• And updated after the first imputation with the imputed data.

• \(\epsilon\) is normally distributed residual error

• Imputed value: \(Y_{imp} = \dot{\beta}_0 + X_{mis}\dot{\beta}_1 + \epsilon\)

• Parameters \(\dot{\beta}_0\) and \(\dot{\beta}_1\) are drawn from their posterior distribution.

• \(\epsilon\) is normally distributed residual error

• Estimate the predicted value using Bayesian regression imputation
  • \(Y_{imp} = \dot{\beta}_0 + X_{mis}\dot{\beta}_1 + \epsilon\)
  • Parameters \(\dot{\beta}_0\) and \(\dot{\beta}_1\) are drawn from their posterior distribution.
  • \(\epsilon\) is normally distributed residual error
• Select \(k\) nearest neighbors to this predicted value from the observed data
• Randomly draw one donor to use as imputed value

• Bayesian logistic regression: logreg
• Bayesian polytomous regression: polyreg
• Classification and regression trees: cart

• Imputed value: \(Log\frac{P(Y_{miss})}{1-P(Y_{mis})} = \dot{\beta}_0 + X_{mis}\dot{\beta}_1 + \epsilon\)

• Parameters \(\dot{\beta}_0\) and \(\dot{\beta}_1\) are drawn from their posterior distribution.

• \(\epsilon\) is normally distributed residual error

• Fit a multinomial regression model.

• Parameters are drawn from their posterior distribution (Bayesian).

• Compute the predicted category.

• Add normally distributed residual error to account for sampling variance.

For continuous or categorical variables

• Draw a random bootstrap sample as training set.
• Fit the tree model using the observed data.
• Find the predicted terminal node for the missing value
• Use the observed cases at the predicted terminal node as donors
• Randomly draw the imputed value from the observed donor cases

library(mice)
help(mice)

Most important arguments in mice()

• m: number of multiple imputations
• method: imputation method
• predictorMatrix: matrix where value of 1 means that the column variable is used as a predictor for the target in the rows.
• maxit: number of iterations
• seed: offsetting the random number generator

summary(nhanes2)

>     age          bmi          hyp          chl       
>  20-39:12   Min.   :20.40   no  :13   Min.   :113.0  
>  40-59: 7   1st Qu.:22.65   yes : 4   1st Qu.:185.0  
>  60-99: 6   Median :26.75   NA's: 8   Median :187.0  
>             Mean   :26.56             Mean   :191.4  
>             3rd Qu.:28.93             3rd Qu.:212.0  
>             Max.   :35.30             Max.   :284.0  
>             NA's   :9                 NA's   :10

plot(imp)

stripplot(imp)

xyplot(imp, chl ~ bmi)

fit <- with(imp, lm(chl ~ age))

pool(fit)

> Class: mipo    m = 5 
>          term m  estimate     ubar         b        t dfcom       df        riv
> 1 (Intercept) 5 180.61667 132.0368  11.14792 145.4143    22 17.69014 0.10131643
> 2    age40-59 5  21.98333 358.3857  15.34230 376.7964    22 19.03236 0.05137138
> 3    age60-99 5  43.78333 396.1105 151.14792 577.4880    22 10.34211 0.45789625
>       lambda       fmi
> 1 0.09199575 0.1797674
> 2 0.04886131 0.1352014
> 3 0.31408014 0.4169004